function theta_prodmult_solve = C_theta_prodmult_fn(theta_prodmult_guess,inputs_struc,fns_struc)

% This function basically just calls another function with an fsolve to
% figure out the theta and prod_mult which match the wage premia and pcGDP
% across counterfactuals
tol = inputs_struc.tol;

% Function handle which gets inputed into the fsolve.
theta_prodmult_zero_handle = @(theta_prodmult)theta_prodmult_zero(theta_prodmult,inputs_struc,fns_struc);
% theta_prodmult_zero_handle(theta_prodmult_guess);

% Solves for the theta and prod_mult.
options_q = optimset('Display','off');
% options_q = optimset('TolFun',tol*10,'TolX',tol*10);
theta_prodmult_solve = fsolve(theta_prodmult_zero_handle,theta_prodmult_guess,options_q);
% theta_prodmult_solve = fzero(theta_prodmult_zero_handle,theta_prodmult_guess,options_q);
  

end

function error_theta_prodmult = theta_prodmult_zero(theta_prodmult,inputs_struc,fns_struc)

w_U = inputs_struc.w_U; w_S_target = inputs_struc.w_S_target;sigma_us = inputs_struc.sigma_us;
prod_grid_size=inputs_struc.prod_grid_size; no_q = inputs_struc.no_q; A_S = inputs_struc.A_S; 
price_size_slope_target = inputs_struc.price_size_slope_target;
A_q_i = inputs_struc.A_q_i; size_q_target = inputs_struc.size_q_target;
p_q_i_calib = inputs_struc.p_q_i_calib; pdf_A_q_i = inputs_struc.pdf_A_q_i;
N = inputs_struc.N; pcgdp_calib_own_prices = inputs_struc.pcgdp_calib_own_prices;
pcgdp_target = inputs_struc.pcgdp_target; theta_q_i = inputs_struc.theta_q_i;

MC_L_part_fn = fns_struc.MC_L_part_fn;

prod_mult = theta_prodmult(1);

price_ratio_target_q  = exp(price_size_slope_target*log(size_q_target));  
price_ratio_target_q = price_ratio_target_q/price_ratio_target_q(1,1);      % Target for aggregate prices across qualities.
MC_L_part = MC_L_part_fn(w_U,w_S_target,theta_q_i(:,1),sigma_us,A_S);
A_scale_q_bisec = (MC_L_part(1,1)./MC_L_part)./price_ratio_target_q; %ones(no_q,1); %
A_q_i_bisec = repmat(A_q_i(1,:),[no_q,1]).*repmat(A_scale_q_bisec,1,prod_grid_size);
inputs_struc.A_q_i = A_q_i_bisec*prod_mult;
Output_struc = solve_equi(inputs_struc,fns_struc);


d_q_i_ch_theta=Output_struc.d_q_i;
M_q_ch_theta=Output_struc.M_q;
pcgdp_prodvec_ch_theta_calib_prices = sum(sum(p_q_i_calib.*d_q_i_ch_theta.*pdf_A_q_i.*repmat(M_q_ch_theta,[1,size(p_q_i_calib,2)])))/N;
pcgdp_prodvec_ch_theta_rel_calib_prices = pcgdp_prodvec_ch_theta_calib_prices/pcgdp_calib_own_prices;
error_prod_mult = (pcgdp_prodvec_ch_theta_rel_calib_prices/pcgdp_target - 1);
        
error_theta_prodmult = [error_prod_mult];

end
